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1 .  INTRODUCTION 

Our  research  has  been  progressing  in  two  directions: 

1)  A  traveling  wave  approach  directed  to  furthering  the 
results  as  in  e.g.,  [10],  [11],  [12].  In  particular,  we  have 
concentrated  on  the  problem  of  identification  of  dominant  modes  of 
vibrations  (dominant  traveling  wave  frequencies)  [8],  In  [8],  we 
have  established  an  analogy  between  vibrations  traveling  along  a 
flexible  structure  and  signal  returns  from  targets  to  an  array  of 
sensors . 

2 

2)  H  /H  theory  of  infinite  dimensional  systems  along  the 
lines  of  the  results  we  had  established  in  [1-3]. 

Our  work  in  both  of  these  is  still  in  progress,  and  will 
continue.  These  are  described  in  more  detail  in  Sections  2  and  3, 
respectively. 
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2. 


CONTROL  OF  FLEXIBLE  LARGE  SPACE  STRUCTURES 


As  it  is  evidenced  in  many  recent  publications  on  the  subject 
(e.g.,  many  papers  in  several  sessions  of  recent  conferences),  the 
main  theme  of  research  is  on  Identification  and  Control  of 
finitely  many  dominant  modes  of  vibrations,  using  finite  dimen¬ 
sional  compensators.  Whether  this  be  in  State-Space  form  or  in 
Frequency  Domain  (e.g.,  considering  those  structures  whose 
transfer  matrices  fall  into  the  Callier-Desoer  Algebra) ,  it 
requires  the  modeling  of  a  finite  number  of  the  (infinitely  many) 
modes  of  vibrations,  and  thus  involves  an  approximation  of 
Hyperbolic  Partial  Differential  equations  (see,  e.g.,  [10]).  The 
modes  to  be  modeled  for  control  must  satisfy  the  differential 
equation  describing  the  structure.  A  common  procedure  for  this  is 
to  consider  hyperbolic  differential  equations  of  the  form 


(1) 


2 

a  u(t,x) 


2 

at 


au(t,x) 

+  A  - 

1  at 


+  A  (x)u(t,x) 
2 


f (X, t) 


where  A^(x)  is  a  partial  differential  operator  in  x  (usually  the 
second  or  the  fourth  order) .  Using  an  orthonormal  set  of 
functions  4>  (x)  which  are  the  eigenfunctions  of  A^(x),  (1)  can 
be  decoupled  into  an  infinite  set  of  ordinary  differential 
equations  for  the  (temporal  oscillation)  modes  of  (1) ,  each 
associated  with  one  of  ^ (x)  (that  describes  the  spatial  amplitude 
distribution  of  the  mode)  (See  e.g.,  [4],  [9]).  Among  these 
infinitely  many  modes,  using  a  finite  dimensional  compensator,  we 
can  only  control  a  finite  number  of  ones. 
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These  can  be  then  described  in  the  usual  state-space  form  as 

(2)  X  =  Fx  +  Gu 

(3)  y  =  Hx  +  +  V. 

In  (2)  ,  F  describes  the  dominant  modes  of  vibrations  we  wish 
to  control,  G  describes  the  distribution  of  control  inputs,  and 
u  is  the  control  input.  It  is  important  to  note  here  that  in  our 
approach  to  identification  of  the  dominant  modes,  we  do  not  need 
to  assume  any  a  priori  knowledge  of  the  frequencies  of  vibration. 

In  (3)  y  is  the  total  of  the  measurements  from  the  struc¬ 
ture;  V  is  the  measurement  noise;  H  describes  the  distribution  of 
the  sensors,  and  v^  describes  the  effect  of  the  modes  that  are  not 
modeled  in  F. 

The  problems  that  are  being  considered  are 

i)  identification  of  the  dominant  modes  (i.e.,  the  unknown 
F)  , 

ii)  choosing  the  right  sensor  structure  and  then  identifi¬ 
cation  of  H, 

iii)  choosing  the  right  distribution  of  actuators  (determina¬ 
tion  of  G) 

iv)  finding  the  right  controller  structure  to  provide  enough 
damping  to  the  dominant  vibrations  despite 

a)  the  measurement  noise  v, 

b)  the  unmodeled  modes  v  . 

1 

Particularly  important  is  the  consideration  of  v^,  as  it  has 
been  shown  that  it  can  cause  instabilities  ([4]).  This  is  refer¬ 
red  to  as  Spillover.  The  noise  v  can  be  handled  considerably  by 
usual  stochastic  control  techniques. 
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It  appears  that  certain  system  identification/signal  process¬ 
ing  techniques  have  already  been  studied  for  applicability  to  this 
problem.  For  example,  ARMA  models  and  their  identification  [5,7], 
frequency  filtering  for  spectral  decomposition  via  Discrete  Four¬ 
ier  Transform  [5],  and  for  isolation  of  the  control  input  from  v^ 
(to  prevent  spillover)  [4].  It  is  recognized  that  [5],  the  domi¬ 
nant  modes  must  be  identified  on-line.  Certain  such  techniques 
have  already  been  in  use  in  approximation  of  partial  differential 
equations  for  a  long  time  [10]. 

We  have  recently  concentrated  on  System  Identification, 
Frequency  Filtering  and  finding  the  right  location  of  the  sensors 
for  this  purpose.  We  have  established  (see  [8],  Appendix  A)  an 
analogy  with  arrays  of  sensors  used  for  target  signal  returns 
(especially  those  for  Radars) .  As  a  consequence,  we  have 
developed  several  new  results  in  ([8],  Appendix  A).  It  turns  out 
that  suitable  location  of  sensors  can  be  helpful  simultaneously 
for  identification  of  the  dominant  modes  as  well  as  in  isolation 
(filtering)  of  the  unmodeled  modes  from  the  output  to  prevent  the 
Spillover.  These  are  described  in  some  detail  in  [8]  (Appendix  A, 
Section  5) . 

i)  The  system  identification  problem  we  have  considered  in 
[8]  (Appendix  A)  is  for  systems  with  zero  inputs.  We  have  been 
working  on  extending  these  techniques  to  systems  with  possibly 
nonzero  inputs,  and  to  time-varying  systems. 

ii)  We  have  been  investigating  the  Toeplitz  (and  related) 
matrices  that  arise  in  the  covariances  in  our  technique.  This 
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would  enable  us  to  approximate  the  dominant  modes  among  infinitely 
many  modes  of  vibrations  in  a  more  effective  way. 

iii)  We  have  been  investigating  the  dualization  of  the  results 
on  sensor  location  to  actuator  locations  for  use  in  disturbance 
cancelling,  and  extension  of  our  system  identification  techniques 
to  control  of  flexible  structures. 

Our  approach  here  has  been  traveling  wave  modeling  of  the 
vibrations,  and  the  analogy  we  have  established  with  arrays  of 
sensors  arising  in  target  signal  returns,  taking  into  consider¬ 
ation  in  full  the  structure  of  the  partial  differential  equations 
that  arise.  This,  in  fact,  appears  to  be  a  very  general  approach 
to  systems  described  via  partial  differential  equations  (see, 
e.g.,  [12]). 
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2  ® 

3.  H  /H  THEORY  OF  INFINITE  DIMENSIONAL  SYSTEMS 


The  problem  we  have  been  working  on  is  the  following. 

Let  G  be  a  semi-infinite  matrix  whose  columns  are  in  i  . 

2 

Let  F  be  linear  operator  on  i  .  Consider  the  sequence 


2  k  ® 

{G,  FG,  F  G . )  =  {F  G) 

k=0 


Now  consider  a  generalized  Gram-Schmidt  process  whereby  we 


consider  the  orthogonal  direct  sums 


k-1 

X  =  ImG  +  ImFG  +  . . .  +  ImF  G 
k 


=  Y  ©  Y  ©....©  Y 

0  1  k-1 


where  the  subspaces  Y  are  mutually  orthogonal.  This  is  similar 

j 

to  the  idea  of  a  Lancoz  Sequence  where  for  a  finite  vector  x  and  a 
finite  matrix  A,  the  sequence 


2 

{X,  Ax,  Ax,  . . . . } 


is  orthogonal ized  via  Gram-Schmidt.  In  this  case,  the  orthogonal 

k  ’  k-1 

projection  of  A  to  the  space  spanned  by  {x,  Ax,....,  A  ) 

X  X 


can  be  recursively  obtained.  This  is  related  to  orthogonal  poly¬ 
nomials  (Lancoz  polynomials) . 

A  similar  technique  has  been  obtained  in  [7]  in  the  incase 
when  F  is  the  shift  operator. 


We  are  working  on  extension  of  this  to  operators  (F,G)  where 
(F,G,H)  is  the  realization  [1]  of  a  transfer  matrix 
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-1  -1 
H  =  Q  R  =  P  Q 
1  1 


where  Q,  R, 


P  Q  are  H  matrices  Q,  Q 
11  1 


being  admissible 


[1]  . 


Such  an  extension  applied  to  the  reachability  map  [G,FG,....] 
yields  a  finite  dimensional  approximation  to  H.  One  can  recur¬ 
sively  obtain  higher  order  approximations.  This  would  also  have 
application  to  Least  squares  identification  of  a  linear  system, 
and  is  also  related  to  Lattice  filters.  This  research  is 
presently  continuing. 
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ABSTRACT 


In  this  paper,  we  give  a  new  extended  an'^  unifying  system 
identification  technique  for  a  class  of  systems  that  include  all 
main  signal  models  that  arise  in  the  Hannonic  Decomposition 
Problem.  This  technique  unifies  and  extends  the  previously 
developed  system  identification  techniques  which  are  improvements 
on  the  PISARENKO  (MUSIC,  dually)  Harmonic  Decomposition  as  they 
arise  in  arrays  of  sensors.  The  advantages  of  the  technique  and 
some  of  its  specializations  given  here  include  having  no 
assumptions  of  stationarity  on  the  stochastic  processes  involved. 
Another  contribution  of  this  technique  is  to  show  that  it  can  also 
be  used  without  any  resort  to  probability  theoretic  concepts,  thus 
bypassing  the  approximation  of  autocorrelations  via  time  averages, 
yielding  the  system  parameters  exactly.  This  technique  can 
be  utilized  to  determine  the  dominant  modes  of  vibrations  of 
flexible  structures  as  well.  An  analogy  is  established  between 
arrays  of  sensors  for  target  signal  returns  and  those  that  can  be 
used  for  vibrations  in  flexible  structures. 
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1.  INTRODUCTION 


In  this  paper,  we  establish  a  new  extended  and  unifying  sys¬ 
tem  identification  technique  for  a  class  of  systems  that  include 
all  basic  signal  models  that  arise  in  the  Harmonic  Decomposition 
Problem. 

The  problem  considered  is  the  following: 

Given  the  system  model 

n 

X  =  Fx  ,  ;  x(0)  e  R  , 

t+1  t 

(1) 

y  =Hx  +v,  t=0,  1,  2,... 

t  t  t 


and  the  output  sequence  y  y  ;  determine  H  and  F. 

n  p  ,  , 

Here  x  e  R  ,  and  y  e  R  .  v  is  the  measurement  noise,  with 
t  t  t 

E  =  R^.  The  only  assumption  on  that  it  should 

decorrelate  after  a  finite  time  a  >  1  (this  assumption  is  also 
used  in  [3]).  This  signal  model  includes  most  of  the  signals 
arising  in  the  Harmonic  Decomposition  Problem  (e.g.,  determining 
the  frequency  spectrum  of  a  given  signal,  as  in  arrays  of  sensors, 
or  as  in  vibrations  arising  in  certain  flexible  stiructures 
possibly  with  some  damping) . 

In  Sections  2  and  3,  we  develop  our  general  system 
identification  technique.  In  Sections  4  and  5,  we  describe  the 
applications  to  two  important  situations:  i)  Arrays  of  Sensors, 
ii)  Vibrations  arising  in  Flexible  Space  Structures.  Also,  we 
compare  our  technique  to  previously  developed  techniques  such  as 
MUSIC,  PISARENKO  HARMONIC  DECOMPOSITION,  and  their  previous 
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extensions  that  can  be  found  in  the  literature  (See  e.g.,  [2,3] 
and  the  references  there) .  The  advantages  of  the  technique  and 
some  of  its  specializations  given  here  include  having  no 
assumptions  of  stationarity  on  the  stochastic  processes  involved. 
Another  contribution  of  this  technique  is  to  show  that  it  can  also 
be  used  without  any  resort  to  probability  theoretic  concepts, 
bypassing  the  approximation  of  autocorrelations  via  time  averages. 
Thus  it  yields  the  system  parameters  exactly. 

2.  IDENTIFICATION  FROM  OUTPUT  COVARIANCES  (PROBABILISTIC  APPROACH) 
Consider 


(2) 


and 


t-N-a 


t-<7 


t  >  N  +  a  , 


corresponding  to  the  signal  model  in  (1). 
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Using  (1),  (2)  can  be  written  as 


(3) 


and 


(4) 


—  — 

—  — 

V 

H 

t 

t 

HF 

V 

= 

X  + 

t+1 

• 

• 

t 

• 

N 

t+N 

HF 

V 

t+N 

_ 

L,  J 

t-N-«7 


t-N-cj 


t-<7 


Vt 


H 

HF 

• 

• 

X  + 

• 

t-N-a 

N 

HF 

t-N-a 


t-cr 

V 

t-N- 


Assuming  a  >  1  is  chosen  large  enough  so  that 

are  uncorre] 
a 

process,  we  obtain 


V  are  uncorrelated,  and  assuming  that  {v  )  is  a 

t-N-C7  t 


(5) 


E{Y  Y  )  = 

t  t-N-a 


H 

HF 


HF 

0 

N 


•N 


T  T 

E{X  X  }n 

t  t-N-a  N 


-13- 


t 


V  and 
t 

zero-mean 


If 


P  =  E{x  X  }  ; 
t  t  t 


f  If***/ 


then. 


E{x  X 

t  t-N-<7 


N+«7 

}  =  F  P 


t-N-<7 


;  t-N-a  >  0. 


Thus, 


T  t  t— N-t7  iji 

E{Y  Y  }  =  n  F  P  (F^)  n  =  K  ;  t-N-a  >  0, 

t  t-N-a  NO  N  N 


It  is  straightforward  to  see  that  the  outputs  produced  by  a 

system  of  the  form 

X  =  Fx  , 
t+1  t 

n 

y=Hx;  xeR,  t>0, 
t  t  0 

satisfies  a  difference  equation  (i.e.,  has  an  Autoregressive  (AR) 


model) 


Ay  +...+Ay  =0;t>r, 

r  t  0  t-r 


for  every  x  e  R  ,  if  and  only  if 
0 

[A  ,...,A  ]fl  =  0, 

0  r  r 

where  A  's  are  pxp  matrices.  Thus  we  obtain 
i 

7.  THEOREM. 

If  F  is  nonsingular  (no  zero  eigenvalues) ,  and  nonsing¬ 

ular,  then  the  AR  models  for  (1)  are  characterized  as  solutions  of 


AK  =  0, 
N 


or,  equivalently. 


T 

K  A  =  0, 
N 
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where 


A  —  fA  A  J» 

0  N 

It  is  clear  that  usually  in  system  identification,  we  are  not 

interested  in  finding  just  any  AR  model,  neither  all  of  them.  All 

the  AR  models  other  then  minimal  ones  will  contain  redundant  modes 

(eigenvalues  of  F  in  (6)),  that  do  not,  in  fact,  exist  in  the  true 

system  that  produces  the  data  y  ,...,y  .  The  same  is  true  for 

0  M 

the  models  (6),  or  (l) .  Assuming  an  irredundant  state-space  model 

in  (1),  (F,H)  is  an  observable  pair.  Then,  it  is  well-known  that 

the  minimal  (irredundant)  AR  models  are  characterized  via  the 

observability  indices  of  (F,H) . 

Contrary  to  the  usual  assumption  in  most  literature,  for 

multi -output  systems,  in  the  minimal  AR  models,  A  is  not  I  .  From 

r  p 

the  results  obtained  above,  and  from  the  elementary  realization 
theory,  it  follows  that  a  minimal  AR  model  can  be  obtained  by 
finding  the  linear  dependence  relationships  among  the  rows  or  the 
columns  of  the  matrix  K^.  r  is  the  largest  observability  index, 
and  N  must  be  chosen  to  be  greater  than  or  equal  to  the  largest 
observability  index.  This  can  be  done  exactly  the  same  way 
that  canonical  forms  for  a  pair  (F,G)  or  (H,F)  are  obtained  (see, 
[1]/  for  a  detailed  exposition  of  this).  This  leads  to 
minimal  AR  representations  of  the  form 
Q(D)y^  =  0;  t  >  0  , 

where  Q(D)  is  a  pxp  polynomial  matrix 

r  r-1 

Q(D)  =AD+A  D  +...+A,  r>N, 
r  r-1  0 
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D  being  the  advance  operator.  In  particular,  a  minimal  Q(D)  can 

be  computed  which  is  row  proper  with  row  degrees  equal  to  the 

observability  indices,  from  which  a  minimal  state-space 

realization  of  the  system  producing  the  data  y  ,  .  •  •  ,  y  can  be 

0  M 

easily  obtained  (see,  e.g.,  [1])  by  inspection  of  A  's.  Clearly, 

i 

N  must  be  chosen  to  be  greater  than  or  equal  to  the  largest 

observability  index.  Once  a  minimal  (F,H)  is  obtained  as 

described  above  (which  is  unique  up  to  a  nonsingular  state-space 

transformation) ,  it  is  straightforward  to  compute  uniquely  up 

to  the  particular  state-space  transformation,  since  0  is  left 

N 

invertible,  and  (by  assumption)  F  has  no  zero  eigenvalues. 

Once  F,  H,  P  are  determined,  the  autocorrelation  function 
0 

of  v^  can  be  easily  obtained  from 


(7b) 


E  { 


]} 


y 

t+CT 


t  T  t  T 
=  0  F  P  (F  )  0  + 

<7  0  <7 


R 

t 


where 


(8) 


R 

t 


V 

t+a 


T 

V  ]  ) 

t+<7 
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yields  the  autocorrelation  function  of  (Note  that  the  first 

term  in  (7b)  is  known) . 

9.  REMARK. 

Note  that  contrary  to  most  system  identification  techniques, 
our  technique,  in  particular,  has  no  stationarity  assumption  on 
the  processes  (Y^)  fact,  in  many  situations,  the 

eigenvalues  of  F  are  exactly  on  the  unit  circle.  This  violates 
the  conditions  for  stationarity  of  (y^),  even  if  (v^)  were  assumed 
to  be  a  stationary  process.  The  technique  given  above  has  no  such 
assumptions. 

10.  REMARK. 

Here  we  did  not  assume  any  special  structure  on 

1  T 

R  :  =  E  (v  V  }  . 

t  t  t 

1 

If  we  assume  that  i)  R^  has  the  special  diagonal  structure  with 

equal  diagonal  entries 
1 

R  :  =  diag  (r  );  r  >0, 
t  1  1  ~ 


and,  ii)  (v^)  is  a  stationary  white-noise  process,  then  we  can 
directly  consider 


(10a)  Z: 


]  } 


t  T  t  T  2 
=  n  FP(F)n  +  R. 
NO  N 
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2 

Then  R  will  be  a  diagonal  matrix  whose  diagonal  entries  will  be 


the  smallest  eigenvalue  of  Z(=  r  ).  From  an  eigenvalue- 


eigenvector  decomposition  of  Z,  we  can  then  compute 


1  t  T  t  T 

K  :  =  Q  F  P  (F  )  n  . 

N  NO  N 


In  this  case,  we  can  obtain  P  ,  F  and  H  from  K  by  the  same 

0  N 

procedure  that  was  applied  to  explained  in  this  section. 

Then,  this  special  case  becomes  a  generalization  of  the 
well-known  PISARENKO  Harmonic  Decomposition  technique  to  the  case 
where  the  output  of  (1)  is  not  necessarily  assumed  to  be  scalar. 


11.  REMARK. 

It  is  clear  that  the  above  technique  also  applies  to  the  case 
Where  (F,H)  n,ay  have  complex  numbers  as  their  entries,  and  ,v^|, 
{y^}  are  complex  valued.  In  this  case,  we  only  need  to  replace 
"T”  (transpose)  with  (complex  conjugate  transpose)  through¬ 

out,  and  change  the  condition  (7a)  to 

AK  =0, 

N 

or,  equivalently, 

K  A*  =  0  . 

N 

12.  REMARK. 

There  is  a  special  case  where  the  technique  developed  here 

can  be  simplified  considerably.  The  reason  why  y  and  Y 

t  t-N-ff 

are  both  considered  is  that  for  our  technique,  must  have 
linearly  independent  columns.  Suppose  that  the  columns  of  H  are 
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linearly  independent.  Then  consider 


T  t  .  t-a  ^ 

(12a)  E  (Y  y  )  =  n  F  P  (F*)  H*  =:  C  . 

t  t-a  NO  N 


Then,  as  F,  are  nonsingular,  and  since  the  rows  of  H 

are  linearly  independent,  one  can  obtain  (H,F)  and  the  minimal  AR 

representations  of  (1)  via  the  left  kernel  of  C  ,  as  explained  in 

N 

this  section.  Next,  P  can  be  obtained  from  C  .  In  certain 

0  N 

cases,  one  is  interested  in  finding  a  nonzero  row  vector  a  such 
that 

ay^  =0;  t=0,  1,  ..., 

for  all  initial  conditions.  (See  section  4) .  Clearly,  this 
corresponds  to  the  condition 
aH  =  0. 

Obviously,  such  a  vector  a  can  be  obtained  by  computing  a 

nonzero  vector  a*  which  is  in  the  right  kernel  of  C  ,  i.e., 

N 

C  a*  =  0. 

N 

If  we  are  only  interested  in  finding  such  a  nulling  vector, 
it  is  enough  to  consider  the  right  or  left  kernel  of 

1  T  t  .t-a  * 

C  :  =  E(y  y  )  =  HF  P  (F*)  H*  . 

N  t  t-a  0 

i.e., 

aH  =  0  , 
if  and  only  if 

1 

aC  =0, 

N 

if  and  only  if 

H*a*  =  0, 
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if  and  only  if 


C^a*  =  0. 
N 


3.  NONPROBABILISTIC  APPROACH 


There  are  standard  techniques  of  estimating  P^  such  as 


1  M  t  T  T  t 

P  =  -  y  /3  F  X  X  (F  ) 

0  a  ^  t  0  0 


where  ^^is  a  weighting  factor,  and  a  is  some  constant.  It  is 

clear  that  nonsingularity  of  P  is  equivalent  to  the  reachability 

°  2t 

of  the  pair  (F7,  x  )  provided  ^  =  7  for  all  t,  where  7  is  a 

0  t 

nonzero  number,  and  M  is  larger  than  the  order  of  F. 

It  is  important  to  note  that,  as  obvious  from  the  results  in 
Section  2,  provided  that  P^  in  (13)  is  nonsingular,  (F,H)  and  P^ 
can  be  obtained  without  having  to  use  any  probability  theory 
concepts.  In  this  latter  case,  the  autocorrelations  and 
covariances  of  all  the  quantities  are  their  time-averages  as  in 
(13) .  Thus,  then  the  assumption  on  {v^}  becomes  that  its  time- 
average  autocorrelation  function  becomes  exactly  zero  after  some 
finite  time.  Thus,  everything  described  so  far  can  be  carried  out 
(with  no  approximation,  as  opposed  to  the  Probabilistic  Approach) 
in  terms  of  time-averages.  Thus,  the  technique  developed  here  is 
perfectly  applicable  to  many  realistic  situations  where  there  is 


usually  no  data  available  other  than  y  ,  y  ,  .  •  .  ,  y  . 

0  1  M 
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4.  HARMONIC  DECOMPOSITION  FOR  ARRAYS  OF  SENSORS 


There  is  a  substantial  amount  of  literature  on  Harmonic 

Decomposition  oriented  system  identification  techniques  (see, 

e.g.,  [2-3]  and  the  references  there).  Among  these,  [2-3] 

appear  to  us  as  among  the  representative  generalizations  of  the 

PISARENKO  Harmonic  Decomposition  technique. 

In  this  case,  usually  target  returns  are  narrowband  signals 

with  center  frequencies  w  ,  .  .  .  ,  w  .  Then  F  has  the  form 

1  n 

jo)  .At 

F  =  diag  {e  i  }, 

where  At  is  the  sampling  interval.  The  matrix  H  in  (1)  now 
corresponds  to  the  measurements  by  an  array  of  sensors.  It  has 
the  structure 

1 
w 

n 

• 

P-1 
w 

n 


r  is  the  delay  due  to  spacing  of  the  sensors.  We  have  assumed 
i 

in  (14)  that  the  sensors  are  equally  spaced.  Clearly,  our 

techniques  in  Section  3  are  applicable  even  if  the  sensors  were 

not  placed  with  equal  distances.  The  problem  is  to  determine 

both  /  •  •  •  0  (j)  },  and  (r  ,  .  .  .  ,T  ).  The  first  set 
Inin 


(14) 


H  = 


1 

w 


w 


P-1 


1 

w 


w 


P-1 


where 


w 
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yields  the  Doppler  frequencies,  r  's  yield  the  angles  of  arri- 

i 

vals.  In  this  case  the  technique  given  in  Remark  10  coincides 
with  the  technique  developed  in  [2]  for  the  special  case  that 
is  assumed  to  be  a  stationary  white-noise  process.  Clearly,  the 
primary  technique  developed  in  Sections  2  and  3  does  not  have  this 
assumption  on 

A  main  difference  between  the  technique  in  Remark  10,  and 
the  one  given  in  [2]  is  that,  in  [2]  the  target  return  signals  are 
assumed  to  be  stationary  zero-mean  stochastic  processes  with  a 
positive  definite  covariance  matrix.  However,  in  many  cases,  the 
return  signals  are,  in  fact,  deterministic  signals  (sinusoidal 
signals  in  the  case  of  Radars,  for  example) ,  with  the  only 
uncertainty  being  at  the  initial  state  x^.  We  exploit  this  fact 
in  Remark  10.  Thus,  the  technique  in  Remark  10  is  a  directly 
applicable  exact  technique  that  can  be  used  purely  in  terms  of 
exact  time  averages  as  opposed  to  the  technique  in  [2]  which  has 
to  approximate  probabilistic  covariances  via  time  averages, 
introducing  an  additional  approximation  error.  Also,  we  can 
incorporate  the  uncertainty  in  signal  magnitudes  into  x^  since 
F  is  diagonal.  Thus  the  H  matrix  in  (14)  has  a  simpler 
structure . 

Another  difference  between  the  techniques  in  [2]  and  in 
Remark  10  is  that,  in  our  case,  the  determination  of  the  number  of 
targets  is  a  problem  of  determining  linear  dependencies  among  the 
rows  (or  columns)  of  a  positive  semidefinite  matrix,  whereas  this 
is  approached  in  [2]  as  a  hypothesis  testing  problem  via  a 


1 ikel ihood-rat io . 
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In  the  technique  developed  in  Section  2,  the  observability 

of  the  pair  (H,F)  (linear  independence  of  the  columns  of  n  )  is 

N 

possible  via  different  r  *s  for  sources  with  the  same  center 

i 

frequency,  and  these  cannot  exceed  p.  But  the  technique  in 
Section  2  (as  well  as  the  one  in  [2])  are  applicable  to  the 
situation  of  an  arbitrary  number  of  targets  that  can  be  larger 
than  p  (  =  number  of  sensors) .  For  example,  then,  two  sensors 
would  be  enough  for  any  number  of  targets. 

Next,  suppose  that  each  target  return  signal  has  a  different 
frequency,  and  that  p  (  =  number  of  sensors)  is  greater  than 
n  =  number  of  return  signals.  In  this  case  the  columns  of  H  in 
(14)  will  be  linearly  independent.  Then,  the  simplified  technique 


in  Remark  (12)  will  yield  (w  ,  .  .  .  ,w  )  and  {r  ,  ,r  } 

1  n  1  n 

in  terms  of  the  left  and  the  right  kernels  of  C  as  in  Remark  12. 

N 

If  all  the  incoming  signals  to  the  array  are  interference, 


and  if  we  want  to  find  a  row  vector  a  such  that 


ay^  =0;  t=0,l,  ..., 

then  any  nonzero  column  vector  in  the  right  kernel  of  C  will 

N 

do  as  a*,  as  explained  in  Remark  12.  Thus  we  can  both  null  and 

estimate  interference  sources  using  the  left  and  right  kernels  of 

C  .  If  we  do  not  want  to  estimate  w  's  and  r  's  it  is  enough  to 
N  i  i 

consider  only  the  right  kernel  of  C^.  In  fact,  if  we  require 

only  output  nulling,  clearly  we  can  only  consider  C  as  explained 

N 

in  Remark  12.  These  latter  cases  constitute  generalizations  of 
the  techniques  given  in  [3]  even  in  this  case  as  follows.  Suppose 
(H,F)  is  observable,  but  the  columns  of  H  are  not  linearly 
independent.  In  this  case,  the  row  vectors  a  such  that 
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aH  =  0  , 

still  coincide  with  the  row  vectors  a  such  that 

C  a*  =  0. 

N 

Thus,  our  technique  still  works  even  if  the  nuinber  of  interference 
signals  are  greater  than  or  equal  to  the  nuinber  of  sensors  as  far 
as  (H,F)  remains  observable,  contrary  to  the  situation  assumed  in 
[3].  We  should  also  note  that  the  technique  in  [3]  is 
probabilistic,  and  as  such  it  introduces  an  approximation  error  in 
estimating  the  probabilistic  covariances  via  time  averages, 
whereas  our  techniques  can  be  used  directly  via  time-averages  as 
well  if  we  define  as  in  (14),  and  they  do  not  require  this 
extra  approximation,  (see  Section  3).  Also,  note  that  via  our 
development,  the  results  obtained  are  still  applicable  even  if 
sensors  are  not  equally  spaced.  In  the  case  of  equally  spaced 
sensors,  the  output  (interference)  nulling  vectors  a's  have  also 
the  following  interpretation  that  has  been  missed  in  [3].  If 


(15)  a  -  [a  , - ,a  ] 

0  p-1 

is  a  vector  such  that 


(16) 


*  1  * 
Ca=0,  or,  Ca=0, 
N  N 


then,  clearly 

(17)  a  (w  )  =  a  +  a  W  + 

i  0  1  i 


p-1 

+  a  W  =  0. 
p-1  i 
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Thus,  if  a (w)  denotes  the  polynomial  obtained  from  a  vector  a  as 
in  (15-17),  the  roots  of  a  greatest  common  divisor  of  all  such 
polynomials  are  w  ,  .  .  .  ,  w 

1  p-1 

17a.  REMARK 

A  technique  for  finding  this  polynomial  whose  roots  are 
w  w  is  as  follows: 

1  p-1 

Let  {7*,  7*,  •  •  .  ,7*)  be  a  basis  for  the  kernel  of  C 
1  2  r  N 

1 

or  (equivalently,  of  C^) .  Then  all  vectors  a  in  this  subspace 
are  given  by 


I  a  7*, 

k=l  ^  ^ 


for  some  arbitrary  a^,  . 


€  R  (or  C) . 


common  divisor  of  the  polynomials  7^(w),  .  • 
obtained  as  described  above  from  7^,  .  •  .  , 
polynomial . 


Thus  a  greatest 

•  f  7  (w) 

r 

7  is  the  desired 
r 


18.  REMARK. 


It  is  well-known  that  MUSIC  is  the  spatial  dual  of  the 
PISARENKO  technique.  Thus  our  results  can  be  similarly  dualized 
to  obtain  spatial  generalizations  of  MUSIC. 
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5.  VIBRATIONS  IN  FLEXIBLE  STRUCTURES 


In  what  follows,  first  we  establish  an  analogy  between  an 
array  of  sensors  receiving  target  return  signals  after  delays,  and 
the  vibrations  occurring  in  a  flexible  structure.  For  example, 
consider  two  types  of  flexible  structures; 


(19) 


2 

a  u(t,x) 


2 

at 


2 

2  d  u(t,x) 

c  -  (elastic  string) , 

2 

ax 


or, 


(20) 


2 

a  v(t,x) 


2 

at 


4 

2  a  v(t,x) 

c  -  (rod) . 

4 

ax 


We  can  approach  the  rod  problem  (20)  first  considering 


2 

a  v(t,x) 
u(t,x);  =  - 

2 

ax 


This  transforms  the  situation  (20)  to  the  wave  equation  situation 
(19).  Thus  we  will  concentrate  on  the  wave  equation  (19). 

It  is  well  known  that  the  general  solution  of  (19)  can  be 
written  as 

u(t,x)  =  f(t  -  -  )  +  g(t  +  -  ) 
c  c 
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which  represents  two  waves  f(t)  and  g(t)  traveling  in  opposite 
directions  along  the  string  with  speed  c.  The  functions  f  and  g 
are  determined  via  the  boundary  and  the  initial  conditions.  For 
simplicity,  let  us  consider  a  semi-infinite  string.  In  this  case 

u(t,x)  =  f(t  -  -  ) 
c 

which  is  a  wave  f(t)  traveling  in  one  direction.  The  case  where  g 
also  exists  can  also  be  approached  similarly.  But  for  simplicity 
of  exposition  of  the  ideas  here,  we  assume  g  =  0.  A  more  general 
traveling  wave  approach  for  more  complicated  structures  can  be 
found,  e.g. ,  in  [5] . 

The  function  f(t)  which  represents  the  oscillations  at  x  =  0 
is  a  sum  of  sinusoids  (infinitely  many) 

"  jwi,t 
f(t)  =  I  A  e  ^  . 
k=l  k 

Now  we  have  a  situation  quite  analogous  to  the  case  of  an  array  of 
sensors  receiving  f(t)  as  a  signal  emitted  by  targets  after  a 
certain  delay  each.  This  analogy  can  now  be  exploited  to  place 
the  sensors  on  a  vibrating  structure  in  certain  ways  that  makes 
the  problem  of  estimating  the  vibration  modes  much  easier. 

To  explain  this,  first  consider  a  single  sinusoid 
jwt 

s(t)  =  Ae  ;  t  e  R. 

Suppose  that  we  place  sensors  at  equally  spaced  points  along 
the  string,  where  d  >  0  is  the  distance  between  the  sensors 
(although  our  techniques  in  Sections  2  and  3  are  still  applicable 
for  arbitrary  number  of  sensors  spaced  arbitrarily) . 
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Then  each  sensor  measures 

d 

jw(t-r-) 

y  (t)  =  Ae  ^  +  V  (t)  ;  t  e  R;  0  <  r  <  p, 

r  r 

where  v  (t)  is  the  measurement  noise, 
r 

Thus ,  we  have 


(21) 


y^(kAt) 

Y  (kAt) 
P 


=  A 


(j> 


(j) 


•j(p-l)r  u 

e 


jwkAt 

e  +  v(kAt) 


where 


T  #  )C““  0/  1^  2^...^ 

c 


and  At  is  the  sampling  interval.  Suppose  that  there  are  n  domi¬ 
nant  modes  w  's  that  we  want  to  identify,  and  the  unmodeled  modes 
i 

contribute  negligibly.  Then  we  have 

X  =  FX  , 
t+1  t 


where  F  describes  the  dominant  modes  of  vibrations,  and 


1 

w 

1 


y 


t 


X  +  V 

t  t 


P-1  p-1 

w  w 

1  2 


p-1 

w 

n 


noise 


=  Hx  +  v,*  t  =  0,  1,2 
t  t 


9 
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where 


i 


V  =  v(kAt) ;  k  =  0,  1,  2,...  . 

k 


Thus,  now  we  have  a  simpler  system  identification  problem  where  H 

has  a  nice  structure  as  in  Section  4.  Clearly,  F,  H  can  be  found 

via  the  system  identification  techniques  developed  here.  In  fact, 

d 

in  this  case  r  =  —  is  known.  Thus  if  we  place  more  sensors  than 

c 

the  number  of  dominant  modes,  the  much  simpler  technique  given  in 

Remark  12  (specialized  in  Section  4)  can  be  utilized  to  estimate 

1 

w  's,  via  kernel  of  C  ,  or  C  . 
i  N  N 

The  case  where  there  is  damping  of  the  modes  of  vibrations 
can  be  very  similarly  solved  via  these  techniques.  One  can  pursue 
the  analogy  we  have  established  further.  There  are  alternatives 
to  placing  p  different  sensors  as  described  above.  One  can 
consider  placing  a  moving  sensor  along  the  x-axis,  and  effectively 
achieve  the  same  result.  This  is  similar  to  the  idea  of  a 


Synthetic  Aperture  Radar  (SAR) .  Note  that  in  the  case  of  flexible 
structures,  the  problem  is  simpler  because  no  angle  of  arrival 
estimation  is  needed. 


22.  REMARK 


The  analogy  with  array  of  sensors  we  described  above  leads  to 
some  other  results  for  flexible  structures.  For  example,  suppose 
we  weigh  the  output  of  each  sensor  via  a^  and  add  them.  Then  we 
obtain  the  response  to  an  incoming  wave  of  frequency  w  as 


d 

N-1  _jj^— w  jwt  _ 

=  A- (  X  a  e  c  )  e  +  v(t) 
k=0  ^ 


jut 

=  A-H(w)e  +  H(«)v(t)  . 


Suppose,  for  example  that,  to  prevent  spillover  (see  e.g., 

[4]),  we  choose  a  's  in  such  a  way  that  H(ci>)  is  a  low  pass 

k 

finite  impulse  response  filter  whose  passband  includes  the  highest 
dominant  frequency  of  vibration. 

Thus,  the  spillover  problem  [4]  due  to  high  frequency 
unmodeled  modes  would  be  eliminated. 
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